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1 Preliminaries 

Considering an ellipsoid whose principal semi-axes a, b, c are ordered in such 
a way that a ^ b ^ c ^ 0, we take the (b, c)-plane as the 'equatorial' plane 
and axis a as the 'polar' axis. A generic point on the surface of the ellipsoid 
is then represented in spherical coordinates R, O, <J> where R is the 'radius' 
(i.e. the length of the vector positioning the point from the centre), is the 
'co-latitude' and <3> is the 'longitude'. 

Therefore, is the central angle positioning the generic point within the 
meridian ellipse passing through the point, whereas $ is the central angle 
positioning the equatorial projection of the generic point within the equato- 
rial ellipse. Obviously, R will depend on the values chosen for O and 

The central cartesian coordinates of the generic point are: 

X = i?cos© , 

Y = it! sin cos $ , 

Z = R sin sin . 

When one uses the eccentric anomalies 9 and cp for the ellipses mentioned 
above, instead of their central angle counterparts, the coordinates become: 

X = a cos 9 , 

Y = b sin 8 cos <p , 

Z = csin# simp . 

The radius to a point upon the ellipsoid surface is then 



R = VX 2 +Y 2 + Z 2 = J a 2 cos 2 6 + b 2 sin 2 9 cos 2 <p + c 2 sin 2 9 sin 2 ip . (1) 



*formerly: European Space Agency, ESTEC, Noordwijk, The Netherlands. 



One may verify the equation of the ellipsoid in cartesian coordinates: 

x 2 y 2 ^ 

a 2 b 2 c 2 

Expressing the latter equation in terms of the central angles and provides 
another expression for the radius at the generic point: 

abc 



R 



\Jb 2 c 2 cos 2 + c 2 a 2 sin 2 G cos 2 <3> + a 2 b 2 sin 2 6 sin 2 $ 



(2) 



For each point upon the ellipsoid, there is another interesting quantity, 
namely the shortest distance H (call it 'height') from the centre to the tan- 
gent plane to the ellipsoid at the point considered. Hence H = Rcos v, where 
v is the acute angle between the radius and the normal to the ellipsoid sur- 
face at the point. Differentiating the equation of the ellipsoid in cartesian 
coordinates, one obtains 4^ dX + p dY + ^ dZ = 0. This shows that the 
vector with cartesian components [4^, X, 4^1 is orthogonal to the tangent 
plane and therefore normal to the surface; let N be its norm. The equation 
for the ellipsoid also shows that the scalar product of this normal vector by 
the radius vector of cartesian components [X, Y, Z\ must be equal to one; in 
other words, NRcos v = 1. Hence one must have 



X 2 i Y 2 , Z 2 



abc 



(3) 



\Jb 2 c 2 cos 2 9 + c 2 a 2 sin 2 # cos 2 ip + a 2 b 2 sin 2 9 sin 2 9J 
abc\Jb 2 c 2 cos 2 + c 2 a 2 sin 2 cos 2 3? + a 2 b 2 sin 2 @ sin 2< 3? 
\/6 4 c 4 cos 2 8 + c 4 a 4 sin 2 9 cos 2 <£> + a 4 6 4 sin 2 6 sin 2 $ 



2 Surface Integral 



Taking here the eccentric anomalies as variables, the general expression for 
the value of a surface integral is (see e.g. Ref.[lJ and Section [6j): 

S= / ^U{9, if)V{9, ip) - W 2 (6, ip) d6 dip where 
Jo Jo 



U{9,<p) 

V(9,<p) 
W(9,<p) 



fdX\ 2 fdY\ 2 (dz\ 

Km) + {m ] + ( ~ ] 



dX\ 2 fdY^ 1 
dip ) \d<p> 



89 J 



dZ\ 

' 1 W 



dX dX dY dY dZ dZ 
d9 dcp d9 dtp d9 dip 
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Performing the calculation gives 



n7T I 
\/b 2 c 2 cos 2 9+c 2 a 2 sin 2 9 cos 2 ip+a 2 b 2 sin 2 6 sin 2 <p s'm9 d9 dip. (5) 
v 

One can perform a cyclic permutation of the semi-axes a, b, c in Equ.([5]). This 
comes down to choosing the (c, a)-plane or the (a, 6)-plane as the equatorial 
plane; the angles 9 and (p belong of course to different planes in each case. 
We found that the easiest way to calculate the surface integral necessitates 
the choice we have made; this, to avoid introducing complex numbers and, 
later on, elliptic functions with parameter greater than one or negative. 

On the other hand, for expressing the surface integral, one can also use 
the central angles (0,$) and obtain, after calculation, an alternate form to 
Equ.(|5]), be it of more complicated nature: 

S = a 2 b 2 c 2 * 

2lT r ^b^c A cos 2 + c 4 a 4 sin 2 6 cos 2 $ + a 4 6 4 sin 2 G sin 2 $ „ ,^ 1JK 
[b 2 c 2 cos 2 + c 2 a 2 sin 2 cos 2< 3? + a 2 b 2 sin 2 sin 2 <I>] 

(6) 




o J o 



3 Other Theoretical relationships 

Equ.(|3]) shows immediately that Equ.(|5]) for the area is in fact 

S = abc I I — — — - sin 9 d9 dip 

H{9,ip) 



Jo Jo 

However, in Equ.([5]), 9 and ip are dummy variables; they may very well be 
interpreted as being O and since the latter have the same intervals of 
integration as 9 and ipl By doing so, the surface integral may be written 
under the alternate form, after Equ.([2]), 

S = abcl / — — — sin©dOd$. 



o Jo 



i?(e,$) 



If one notices that the volume of the ellipsoid is T = Airabc/3 and that the 
infinitesimal element of solid angle d^ in the direction (O, <I>) is, in spherical 
coordinates, sinOdOd^, one obtains the remarkable result that 

S f 4n 1 df 



T J R(^) 4vr ' 

where ^ is the value between and Air of the solid angle positioning the 
generic point of radius R. 
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The equivalent at two dimensions of the latter expression is the ratio between 
the length A of an ellipse and its area A = irab, using the central radius and 
the central angle: 

A _ f 2n 1 d$ 
A~ 2 J R{¥)^ ' 

Finally, consideration of Equs.([2]),(j3]) and ([6]) shows that one can also write 

S = / u/U W sin9ded$ = / — TTvd*. 
Jo Jo H 9,$ J H * 



4 Calculating the surface integral 

Our starting integral will be Equ.([5]). If one notices that, by symmetry, 
the ellipsoid can be decomposed into eight similar octants, we can limit the 
integration range for 8 and p to [0, ^] and write 



S~ Sab I 2 /" \ i + 



Jo 



a' 



b 2 



+ 1 



b 2 



sm p 



sin 2 9 sin 9 d6 dp 



We first let k(<p) 



these quantities remaining within the interval [0, 1]. We obtain 



S = 8ab 



2 r 2 



o Jo 



I — e 2 + k 2 (p) sin 2 9 sin d0 d</3 



In order to explicit the inner integral, we substitute £ = k cos 9 and get 



Vl-e 2 + k 2 sm 2 9 sm9d9 = y I J (1 - e 2 + k 2 ) - £ 2 d^ . 
o * Jo 



A primitive integral can be found in Ref.[2], §2.262.1: 



1 - e 2 + A: 2 



asm 



V / (l-e 2 + fc 2 ) 



so that 



f, 2 , /2 . 2fl . fl , fl v 7 !^? , 1 - e 2 + k 2 . 

V 1 — e z + Ar sin y sin 9 do = 1 asm 

2 2iJi' 



Vl - e 2 + /fc 2 ' 
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There are two occurrences of k within the argument of the arcsine; for the 
range of values of k and e at hand, we therefore better replace 

k k 
asin , = atan ■ 



Vl - e 2 + k 2 VT^e 2 ' 

The surface integral becomes thus 

S = 2-Kbc + 4ab / . . atan - dw 

Jo H<p) 

2 

or, in terms of variable k, using kdk = (1 — simp cos tpdcp, 



S = 2itbc + Aab / v 7 dk . 

y*- (1+^)1/(1+^) -s 

In order to simplify this form, let r] = 1 + (^r) 2 #s 1, so that 

S = lube + 26c ^atan^FT ^ (?) 



Expressing this integral in terms of known integrals will necessitate an inte- 
gration by parts: 

V 



let u(r/) = atan sj r\ — 1 and = 

First, we obtain easily 1/(77) = 2^7^T ' 

Next, we need a primitive integral for v'(rj); one could take f^ 2 f'(C) d£. How- 

"p- 

a 2 

ever, by substracting from the latter the constant J^ 2 v'(Q d£, one obtains 

62" 

another primitive integral, which is available in Ref.[2] at §3.132.5: 

a 2 

v(v) = ~ I t/(C) dC = "2 [p^F^r?), m) + pE(r(rj),m)] (8) 



with 



asm 



c 2 ry 



6 2 



m 



a 2 (6 2 



6 2 (a 2 



and where F(r, m) and E(r, m) designate the elliptic integrals of the first 
and second kind, respectively. Quantities r and m always remain within the 
intervals [0,^] and [0,1], respectively. 
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Result ([8]) can be checked directly by differentiation of v{rf). For that pur- 
pose, observe that 



_ 1( 9F dE 



dij 



1 



and that — — 
dr) 



1 



dr 
dr] 



With some skill, result ([8]) can also be obtained by using the Mathematica™ 
software. 



We may now express that 



a 2 
71 



uv' drj 



uv 



a 2 
71 



77 



a 2 
71 



u'vdr] . 



77 



One has t(t-) = asinO = and t(j^) = asin 1 = ^ . Careful examination 
of the expression for r using the definitions of r] and k shows indeed that 
r = asin cos <p = | — <p . 



Therefore 



uv 



a 2 
71 



2 atan 



b 2 



1 [p -1 K(m) + pL(mj\ , 



where K(m) and L(m) are the complete elliptic integrals of the first and 
second kind, respectively. 

Coming back to the surface integral, one has 



5 = 2nbc + Abe atan \/ -p- - 1 [p x K(m) + pL(m)] + 26c I 



where 



I 



- 2 p 1 F (t(tj), m) + pE (r(r/), m) 



r]^/r] - 1 



dr] . 



With the knowledge that 
a 2 



b 2 



sin 2 r 



and that 



dr] 



-kdk 



b 2 



) sin (p cos if dip 



b 2 



) cos r sin r dr 



one can express integral I in terms of r, after simplifications, as 



c 
62 



2 [p 1 F(r, m) + pE(r, m)] sin r cos r 
(1 — me 2 sin 2 r)y/ 1 — msin 2 r 



dr . 
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Fortunately, integral / is now the combination of two definite integrals which 
can be found in Ref.[2], the first one at §6.113.2, 



F(r, m) sin r cos r 



o (1 — me 2 sin 2 r) \J\ — m sin 2 r 



dr 



1 



VT^T 2 



me 



atan ( \f\ — m ] K(m) — — F(asine, m) 



and the second one at §6.123, 

E(r, m) sin r cos r 



(1 — me 2 sin 2 r) yl — msin 2 r 
atan I yT" — ;// 



:dr 



me 



< \ 7T . 7r\/l - e 2 (l - Vl - "ie 2 
L(m) h(asm e, mj H 



2e 



Replacing all of these results into the surface integral and simplifying every- 
thing finally gives the surface area of the ellipsoid: 



S = 2vr 



c 2 + 



be 2 



\/a 2 —6 



:F(asin e, m) + b\J a 2 — c 2 E(asin e, m) 



(9) 



with e = a / 1 ?r and m 



a 2 (b 2 -c 2 ) 



b 2 (a 2 -c 2 ) ' 

This formula can be found without demonstration in the add-ons to Ref.[3], 
be it with a severe typographical error in versions earlier than Mathematica™ 5. 



5 Different kinds of ellipsoids 

Interpretation of the parameters e and m is important. Eccentricity e is the 
largest eccentricity of all possible ellipses obtained by cutting the ellipsoid 
by planes. If we let m = sin 2 7, then 7 is the angle between the (b,c)-plane 
and any of the two planes containing the circular sections of the general 
ellipsoid. The circular sections are orthogonal to the (a,c)-ellipse. 

When 7 < atan ^, the ellipsoid shows a prolate character. 

When 7 > atan |, it shows an oblate character. 

When 7 = atan | , the two planes containing the circular sections 
of the ellipsoid cut the (a,c)-ellipse along conjugated diameters 
(see Section [6]) and the ellipsoid shows a spheroidal character. 
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In order to derive all possible limit cases, let us rewrite formula (|9|) as a 
function of the three independent parameters a, e and m only: 



S=2ira 2 



VT 



v 7 ! — me 2 



n o r, n o.F(asine,m) _, 

VI - e 2 V 1 - me 2 +(1 - e 2 )—^ — '- + eE( 



asme, m) 



(10) 

• For m = with < e < 1, one has b = c and the ellipsoid is prolate 
of revolution (about axis a). Since F(asine,0) = E(asine,0) = asine, the 
surface area (|10p simplifies into 







•Sprol.rev. = 27Ta 2 \/l ~ e 2 


/ T asin e 

V 1 - e 2 H 

e 



(11) 



One may verify that 

'S'prol.r 



e^O 



> 4ira 2 (sphere) 



5, 



prol.rev. 



e->l 



> (bar-like shape along a) 



• For m = 1 with < e < 1, one has a = b and the ellipsoid is oblate of 
revolution (about axis c). With the knowledge that J dip = atanhsint/?, 
one has F(asine, 1) = atanhe and E(asine, 1 
(jlOP simplifies into 



cos ip 

e, so that the surface area 



5obl.i 



2na z 



1 + (1 



atanh e 



(12) 



One may verify that 



5oblj 



e->0 



> 47ro 2 (sphere). 



• For both m = 1 and e = 1 though, the factor in front of Equ. (|10p becomes 
indeterminate! But since, by definition, 

VT^ 2 _ b 
Vl — tne 2 a 

this only means that b may be chosen at will in the range [0, a]. Because 
c = 0, we are indeed dealing here with a flat elliptic disc in the (a, 6)-plane. 
Since lime^i [(l — e 2 ) atanhe] = 0, its surface area is 

ffeii.disc = 2-nab , 
the disc becoming circular when b = a. 
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• For e = 1 together with m indeterminate, one must have c = b = 0, so 
that the ellipsoid degenerates into a bar along axis a with, from (I10p . 



Sbar = 



• For e = 0, one has c = b = a and the ellipsoid reduces to a sphere. 
Parameter m is indeterminate in this case but, since 
for all to, one obtains, from (fTU|) . as it should: 



^sphere — 47TQ, 



For completeness, we give in Table [T] a classification of all possible ellipsoids, 
together with their surface area S; we use the notations 

,* la 2 + c 2 , * a 2 

b = \ — - — and m = * 

V 2 a 2 + c 2 



for the intermediate semi-axis and the parameter of the spheroid. 



ellipsoid 


semi-axes 


e 


m 


area S 


prolate of revolution 




< c = b < a 


< e < 1 


TO = 


Equ.QJ 


general prolate 





<c<b<b*<a 


< e < 1 


< to < m* 


Equ-GDU 


general spheroid 





<c<b = b*<a 


< e < 1 


TO = TO* 


Equ.CED 


general oblate 





<c<b*<b<a 


< e < 1 


TO* < TO < 1 




oblate of revolution 




< c < b = a 


< e < 1 


TO = 1 


Equ.CED 


elliptic disc 




= c < b < a 


e = 1 


TO = 1 


27TO& 


circular disc 




= c < b = a 


e = 1 


TO = 1 


27T0 2 


bar 




= c = b < a 


e = 1 


< TO 1 





sphere 




< c = b = a 


e = 


< TO 1 


Ana 2 


point 




= c = b = a 


^ e < 1 


^ TO < 1 






Table 1: Classification of the ellipsoids 



6 Curvature of the ellipsoid 

A description of how to obtain the (two) curvatures at one point of a surface 
is given in Ref.[T]. Using the position vector R in coordinates (9,ip), one 
first needs vectors Rg and R^ which are the derivatives of vector R with 
respect to 6 and <p, respectively. One may then define the (positive definite) 
matrix of metric coefficients 



" u 


w~ 




Re 


w 


V 




Rip 



• [ Re Rip ] ■ 
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0/00 00\ 

After calculation, one finds U = a sin 9 + cos 9 [b cos <p + c sin y?J 

V = sin 2 (ft 2 sin 2 (p + c 2 cos 2 y>) , 
= — (b 2 — (?) sin cos 9 sin </j cos (p. 



Note that the surface area is S = J^TcT II ^ ^ -^vll d#dy? (see Section [2]). 
Next, knowing that the unit vector normal to the ellipsoid surface is 



n 



R e AR 



<p 



RoARu 



one needs its derivatives ng and n^ as well, in order to construct the (sym- 
metric) matrix 







Re 


fJL A 




Rip 



n 



One finds that x 



abc 



A = xsin 
H = 



\/b 2 c 2 cos 2 + c 2 a 2 sin 2 cos 2 <p + a 2 6 2 sin 2 sin 2 </j 

2 /) 



One can now define the curvature tensor 









-1 






' U 


w~ 


Re 




X=[n e n v ] 


w 


V 




Rip 





which is a symmetric tensor because 



(13) 



Re 
Ru 



X- [Re Rip 



fi A 



is a symmetric matrix. 



From ng = x ' Re an d n v = x ■ R v > one concludes that 



dn = X' d-R 



(14) 



Hence, for a given d.R, the more pronounced the curvature, the more vector 
n changes its orientation! 

The (scalar) curvature x f f° r the direction £ = n gj|n is then 



d • x- di? di? • dn 
d.R- di? dR- dR 



[ d9 dp] 


fjL A 




' d9 ' 

dp 


[ d9 dp] 


U W 
W V 


' d9 ' 

dp> 
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When dR is along a principal direction of curvature (eigenvector of x) > then 
dn = x d-R, where the scalar x i s a principal curvature (eigenvalue of x) ■ I n 
such a case, 



R^ 



or equivalently 





' de ' 




Rg 


n e rip] 


dip 


= X 





[Rg R^ 



de 

dip 









' u 


w~ 




' de ' 




"0" 


( 


jJL A 


- X 


w 


V 


) 


dp 








The characteristic equation of this eigenvalue problem allows one to calculate 
the two principal curvatures x 1 an d X 2 • Their sum (twice the so-called mean 
curvature) and their product (square of the so-called Gaussian curvature) 
are found to be, using Equs.(UJ and Q: 



X 1 + X, 



xV + \U - 2fiW H 3 (a 2 + b 2 + c 2 



R 2 



UV -W 2 



i 2 b 2 c 2 



and 



X X c 



xX 



H 4 



UV-W 2 



i 2 b 2 c 2 



(15) 



(16) 



the latter expressions being independent of the choice of coordinates: (0, p) 
or (0,$). Being the eigenvectors of the symmetric tensor the principal 
directions of curvature Z\,£-2 will be orthonormal. 

When, in a general ellipsoid, one translates a circular section parallel to itself, 
one obtains new circular sections whose centres travel in the (a,c)-ellipse 
along a diameter which is said to be conjugated to the diameter through 
which the circular section of radius b passes. When the plane of one of these 
parallel circular sections becomes tangent to the ellipsoid, the circular section 
becomes a single point called an umbilic. There are four such points on a 
general ellipsoid, with cartesian coordinates 



X = ±a 



b 2 



a 2 — c 2 



Y = 



Z = ±c 



b 2 



a 2 — c 2 



so that, for these umbilics, 



R = VX 2 + Y 2 + Z 2 = Va 2 + c 2 -b 2 and H 



X 2 , Y' 2 , Z 2 



ac 

T 



With the use of Equs (|15|) and (|16p , one can calculate the principal curvatures 
at the umbilics: 



ac 



X, 



-1 -2 fc3 ' 

In accordance with the definition given in Ref.pQ, the curvature at an um- 
bilical point is the same in any direction. 
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